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ABSTRACT 

We investigate the interaction of differential rotation and a misaligned mag¬ 
netic field. The incompressible magnetohydrodynamic equations are solved nu¬ 
merically for a free-decay problem. In the kinematic limit, differential rotation 
annihilates the non-axisymmetric field on a timescale proportional to the cube 
root of magnetic Reynolds number ( Rm ), as predicted by Radler. Nonlinearly, 
the outcome depends upon the initial energy in the non-axisymmetric part of 
the field. Sufficiently weak fields approach axisymmetry as in the kinematic 
limit; some differential rotation survives across magnetic surfaces, at least on 
intermediate timescales. Stronger fields enforce uniform rotation and remain 
non-axisymmetric. The initial field strength that divides these two regimes does 
not follow the scaling RrrT 1 ^ 3 predicted by quasi-kinematic arguments, perhaps 
because our Rm is never sufficiently large or because of reconnection. We discuss 
the possible relevance of these results to tidal synchronization and tidal heating 
of close binary stars, particularly double white dwarfs. 


Introduction 


Current understanding of stellar differential rotation leaves much to be desired. 
Centuries of sunspot observations show that the Sun rotates more quickly than its equator. 
Helioseismology reveals that this latitudinal variation extends throughout the convection 
zone, but contrary to expe ctati ons for a (nearly) isentropic region, the angular velocity is 


not constant on cylinders (ISchou et ah 


19981 ). There is no consensus as to how this pattern 


is maintained, though the convection itself is presumably essential. Inversions for the 
radiative core are consistent with uniform rotation down to at least 0.2R G ( Chapl in et ah 
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19991) : constraints at greater depth are weak because they depend on a few low-degree 
p inodes. Asteroseismological analyses of Kepler photometry find that many giants and 
subgiants have cores that rotate more rapidly than their envelop es, but not so ra pidly 


as if these cores had contracted at consta nt angular momentum flDelieuvels et al. 


2012 


Mosser et al 


2012 


likely responsible 


Deheuvels et al 


2014 . 


Maeder fe Mevnet 


Magnetic transport of angular momentum is 


2 014T) . but standard prescriptio ns used by t h e st el 


evolution community fail 


2012 


Cantiello et al. 


o explain the observations quantitatively (jE g gen be rger et al 


ar- 


2014) . White dwarfs, the end-stat es of low-mass star s , are known to 


rotate slowly, with periods ranging from hours to days ((Berger et al. 


2005 


Kawaler 


2014 


and references therein). Asteroseismological attempts to measure differenti al rotat ion in the 


interiors of pulsating white dwarfs have so far produced ambiguous results (ICharpinet et al. 


2009 


Corsico et al. 


2 0111) . The role of stellar magnetic fields in the transport of angular 


momentum in all of these objects may be crucial but is not yet well understood. The 
present work addresses one aspect of this complex problem. 

Ferraro’s isorotation law applied to stellar interior states that “the star can possess 
a steady field only if the field is symmetric about the axis of rotation, and each line of 
force lies wholly in a surface symmetric about the axis and rotat ing with u niform angular 


velocity” (jFerraro 


19371 ). A detailed derivation can be found in ICow ling (119571) . The 


theorem can be mathematically expressed as B p • Vf2 = 0, where B p is the axisymmetric 
meridional field and fl the angular velocity. This formula shows that the contours 
of axisymmetric meridional field are parallel to the contours of angular velocity. The 
isorotation law is strictly valid_only f or a p erfect conductor, i.e. the magnetic Reynolds 


number Rm —> oo. Later, 


Mestel &; Weissl (119871 ) studied the dynamical and resistive 


effects of departures from the isorotation law. Their results suggest that Alfven waves 
transfer angular momentum along meridional field lines—at different speeds on different 
lines—until the isorotation law is achieved. This mechanism of “phase-mixing” was studied 
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by 


Ionson (1978); 


Hevvaerts fe Pr iest 


(1198.4 


Spruit 


( 1999 ). The latter estimated that the 


wave amplitude should decay as e on a phase-mixing timescale t p («7»V 2 ) 1/s , 

R being stellar radius, 77 being magnetic diffusivity, and V K being a typical Alfven speed. 
Therefore, it may be inferred that if all field lines are attached to a uniformly rotating 
solid core, the n the entire stellar (or planetary) in terior will eventually a chieve solid-body 


rotation (e.g., 


Charbonneau fe MacGreg or 


19921) . 


Meste l fe Weissi ( 19871 ) suggested that 


even without such a core, a fluid body will tend toward solid-body rotation if its magnetic 
field is significantly non-axisymmetric: for example, a dipolar field whose axis is not parallel 
to the rotation axis, a situation referred to in our paper as an oblique rotator. It is this 
last assertion, rather than the whole subject of stellar differential rotation, that is the main 
focus of the present paper. 


Mestel and Weiss perhaps overstated their argument by restricting the velocity field to 
pure rotation, rather than a combination of rotation and meridional circulation. Consider, 
for example, a star in a nontrivial state of isorotation aligned with its magnetic axis, but 
undergoing slow precession around another axis (due, perhaps, to the tidal torque of a 
companion). Such a star would not be strictly axisymmetric with respect to the axis defined 
by its total angular momentum, but its velocity field—which would have meridional as well 
as azimuthal components with respect to that axis—could satisfy Ferraro’s Law and be 
steady in a frame rotating with the star’s body axes. The calculations described in this 
paper allow for meridional motions. 


Even as Mestel and Weiss posed the problem, the outcome 


relative strengths of the field and of differential rotation. 


surely 


depends on the 


Radlerl (119861 ) considered the 


kinematic limit in which the backreaction of the magnetic field on the flow is neglected, so 
that the field evolves initially linearly according to the induction equation in a prescribed 
flow. He concluded that the combined effects of differential rotation and magnetic diffusion 
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cause the non-axisymmetric field to decay more rapidly than the axisy mmetri c field, thus 
tending to reduce the magnetic obliquity. Allowing for magnetic forces, 


Meste l (J1999 


9.3) 


concluded that a nonaxisymmetric meridional field with a weak component perpendicular 
to the rotational axis can be destroyed by differential rotation, but a stronger perpendicular 
component can destroy the differential rotation. The principal goal of the present paper is 
to test and quantify these conclusions by explicit calculations. 


The problem posed here has several applications in geophysics and astrophysics. In 
many stars and fluid planets, it will often be complicated by a magnetic dynamo, which 
may reinforce the magnetic field in ways that cannot confidently be predicted. We are 
motivated, however, mainly by applications to white dwarfs, whose magnetic fields are 
probably relics inherited from their progenitors. The interaction between differential 
rotation and magnetic fields may therefore be simpler to study in the context of white 
dwarfs than in main-sequence or red-giant stars. Admittedly the observational constraints 
are weaker in white dwarfs. On the other hand, the consequences might be more dramatic, 
as we now explain. 


In the double-degenerate scenario for Type la supernovae (SNe la), two white dwarfs 
orbiting one another are gradually driven together by gravitational radiation. At late 
phases of the inspiral when the stars are separated by a modest multiple of their radii, they 
will exert mutual tidal torques tending to enforce synchronism between their rotational 
and orbital frequencies. These torques, which likely involve resonant excitation of inertial 
oscillations and internal waves (g modes), are expected to be unequally distributed within 
the star, in fact concentrated toward the surface where thermal time scales are s h ortest . 


densities are 

east. 

Bur kart et al. 

2013 

(Fuller & Lai 

2013) 


east, and tidally excited waves may break nonlinearly ((Fuller fe Lai 


DaH’Osso & Rossi 


2012 ; 


20141). possibly with observable consequences 


20131) . When the star rotates rapidly compared to the tidal frequency, i.e. 
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when hhn in A> — f2 sp m, the torque is also concentrated in latitude, toward the equator 


fjFuller fc Lai 


2 0141) . Magnetic stresses are likely needed to couple the rotation of the stellar 


interior to that of the surface layers. The amount of associated dissipation in the interior 
may depend not only on the strength of the magnetic field and of the tide, but also on the 
symmetry of the field. If it is essentially axisymmetric, then since plasma viscosity is likely 
negligible, some turbulent dissipation is probably required to transport angular momentum 
across the lines. 


Also, the efficiency of magnetic redistribution of angular momentum may affect the 
degree of nonsynchronous rotation, which, even if small, determines the total dissipation 
associated with a given tidal torque. An upper bound to the tidal heating rate is the power 
required to maintain synchronous rotation, E s pin = (A + I 2 )Q(l, where hi is the orbital 
angular velocity and I\^ are the moments of inertia of the two stars. For two 0.7 M 0 carbon- 
oxygen white dwarfs driven together by gravitational radiation, P spin ~ 10 38 P mi ^ 3 ergs^ 1 , 
where P min is the orbital period in minutes, while the time before contact is 400 PA 3 / 3 Y 1 - 
As noted in the works cited above, however, the actual dissipation rate will be less than this 
by an appropriate average of (1 — fl sp i n /hlorb), hl sp in being the rotational angular velocity of 
each mass element (which differs among elements if the star rotates differentially). Thus 
the dissipation will be quite small if the tidal torques are efficiently redistributed so that all 
parts of the star are kept nearly synchronous with the orbit. The dissipation will also be 
small if the tidal torques are weak so that fl spin <C Q or b- Still, if even a small fraction of 
Aspin were dissipated and the heat were transported to the surface, such systems could be 
quite luminous. 


Surface magnetic fields of white dwarfs vary widely. Zeeman me asure ments indicate 


that ~ 10% of these stars have fields exceeding 2 MG flLiebert et ah 


200j3|), but sensitive 


polarimetry (of admittedly small samples) suggest that the majority have surface fields 
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< 10 kG (ILandstreet et al.l 2 01 2). In the expected absence of dynamo action, electrical 
conductivities in the degenerate interi ors o f th ese stars are such that their fields should 
decay on timescales ~ 10 9 yr (jFontaine et al. 19731 ). This is of the same order as the ages 


of most observed white dwarfs and likely also of many SN la progenitors. Depending upon 
the star’s mass, crystallization begins at the center after on e to a few billion years when 


the luminosity has fallen to 10 3 to 10 4 L 0 (IRenedo et al. 


20100- Some field lines may 


be anchored to the growing solid core, an effect not accounted for in the idealized models 
considered here. 

The d oub l e-degenerate scenario for SNe la, though favored by recent evidence, is 


unproven flMaoz et al. 


2 014 ). A variant invokes merging via head-on collis ions in a 


tripl e-star syst em rather than inspiral driven by gravitational waves 


Kushnir et a. 


Brown et al 


20 


3). Yet short-period white-dwarf binaries do exist flNelemans et al 


( Katz fe Dong 


2012 


2004 


20111 ) . and their end states invite speculation even if they are not destined to 


be supernovae. 

The plan of this paper is as follows. §2 frames the free-decay problem and our numerical 


methods. §3 compares results for the kinematic problem, in which the differentia 
is prescribed and the backreaction of the field neglected, with the predictions of 


rotation 


Radler 


019861 ). §4 presents results for the full nonlinear problem, including estimates for the critical 


initial field strength as a function of Rm. Numerical considerations limit our calculations to 
Re = Rm < 10 4 , much less than in a real white dwarf. Finally, §5 discusses the relationship 
of these numerical results to the astrophysical problem that motivates this work. It is clear 
that there are many aspects of the joint tidally-driven evolution of the magnetic field and 
differential rotation that will need further study. 





































2. Formulation 


We calculate the three-dimensional and fully nonlinear unforced MHD equations in 
a spherical shell with inner radius r* and outer radius r a . The differential rotation and 
misaligned field are given as initial conditions. To conserve total angular momentum, 
stress-free conditions are imposed on the fluid velocity and insulating conditions on the 
magnetic held at both boundaries. 


The dimensionless MHD equations for an incompressible fluid of constant and uniform 
density read 

du _ 1 0 . . 

— + u- Vm = -Vp + — -V 2 u + (V x B) x B, (1) 

ot He 

£)T3 1 

— = Vx(uxB) + —V 2 B. (2) 

Length has been normalized to r a , time to Dq 1 , velocity to D 0 r o where Do is characteristic of 
the initial angular velocity [see eq. Q], pressure to pDgr^, and magnetic held to y/]JjIQ 0 r o . 
There are three dimensionless parameters governing the MHD how. The Reynolds number 
Re = where v is kinematic viscosity, is the ratio of the viscous time scale to the 

initial rotational time scale. Similarly, the magnetic Reynolds number Rm = Vt^rl/r], where 
T] is magnetic diffusivity, is the timescale for diffusion of the initial magnetic held relative 
to the rotation time. For numerical feasibility, we keep Re = Rm, namely the magnetic 
Prandtl number Pm is unity. The dimensionless Alfven velocity Va = f?o/(\/PR^o r o), 
where Bq is characteristic of the initial held [see eq. fflOl) ]. measures the strength of initial 
held relative to initial rotation. To minimize the effect of the inner sphere while avoiding 
the coordinate singularity at the origin, we take r^/Vo = 0.1. 


The details of numerical method can be found in 


Hollerbachl ( 120001 ). Toroidal-poloidal 


decompositions are used to guarantee V • u = V • B = 0: 


u = V x (e e r ) + V x V x (/ e r ), B = V x (g e r ) + V x V x {h e r ). 


(3) 
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Note that the toroidal part V x (ee r ) has a latitudinal component proportional to de/dcj), 
and that the poloidal part V x V x (/ e r ) has an azimuthal component proportional 
to d 2 f / d(j)dr . Therefore, when speaking of non-axisymmetric fields, we will refer to 
components parallel to as azimuthal (rather than toroidal), and to components parallel 
to the rO plane as meridional (rather than poloidal), to avoid confusion caused by experience 
with axisymmetry, where “toroidal” is synonymous with “azimuthal,” and “poloidal” with 
“meridional.” 


In spherical coordinates (r, 6 , (ft) the functions {e, /, g, h} are expanded in the angular 
coordinates with spherical harmonics and in radius with Chebyshev polynomials. For 
example, 

e(r, 9 , M = [ e L(P t) Pi m ( cos 9) cos (mcft) + ef m (r, t) P z m (cos 9) sin (me /))], 

l,m 

eL M) = XI (*) T k 0) and e lm M) = e H"» (* ) ^ > ( 4 ) 

A; 

where x = (2r — r 0 — r,)/(r 0 — r,) G [—1, +1]. We use a second order Runge-Kutta scheme 
for time stepping. The diffusive terms are treated implicitly. 


We impose stress-free boundary condition for fluid velocity at both r 0 and r, : , namely 
u r = T r o = T r( p = 0. Translated to spherical harmonics, this becomes 


flm 


d 


1 d \ d 

'Jim I 


C/m ) 0. 


(5) 


dr \r 2 dr Jl " L J dr 

These hold for both cosine and sine components and so the superscripts c and s are omitted. 
We require the field to match onto potential fields interior to r t and exterior to r 0 that are 
regular at r = 0 and r = oo, respectively, so that 

l + l 


9lm 


d 

dr 


him 0 at ?y, Qim 


d l . , 

~r + - ) him = o at r Q . 

dr r 


( 6 ) 


The initial differential rotation profile is taken for the hydrostatic equilibrium, i.e. 
inertial force is balanced by pressure gradient, such that the angular velocity depends only 
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on cylindrical radius R = r sin 6, 


Q = £Iqt 2 sin 2 9. 


(7) 


Nonmagnetic force balance with more general patterns of differential rotation would 
require stratification, which we wish to avoid. It seems unlikely that this simplification 
qualitatively affects the competition between differential rotation and non-axisymmetry, 
but it does restrict the allowable isorotational states in case the former triumphs over the 
latter. A larger exponent of R might better imitate the concentration of tidal torques 
toward the surface and equator but would lead to rapid magnetic and viscous diffusion 
at the numerically accessible values of Re and Rm. Like this profile, the differential 


rotation resulting from tidal torques in an inspiraling binary i s proba 


magnetorotational instability (MRI) because dQ 1 2 /dR > 0 (e.g., 


Balbus 


jlv immune to 


20031). 


For the initial field we choose free-decay modes, eigenfunctions of the induction 


equation with finite conductivity and currents confined to the body (IMoffatt 
the conductivity is uniform, the poloidal expansion coefficients obey 

dhi m d 2 him 1(1 + 1) 


19781 . §2.7). If 


dt 


dr 2 


lj lm • 


( 8 ) 


The eigenfunctions for hi m are linear combinations of spherical Bessel functions of the first 
and second kind that satisfy the insulating boundary condition (|6]). We choose for the 
initial conditions a linear combination of the two lowest-order poloidal free-decay modes 
(l — 1 ,m — 0) and (l — 1 ,m — l)|l 


h(r ) = hio(r) cos a + hu(r ) sin a, 


(9) 


1 These modes suffice to match onto a general dipole field in the vacuum exterior to 

the star. In white dwarfs, modes with higher l and/or more radial nodes would undergo 
significant resistive decay on Gyr timescales. 
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iii which a is the angle between the rotational and magnetic axes. The initial field is then 


B r = 
B e = 

= 


Bq— [h \ 0 cos 9 cos a — sin 9 (hh cos <j) + h s u sin </>) sin a ], 


-Bn- 


B n - 


dhlo . .. „ 

-sm 9 cos a + cos 9 

dr 

dh c n 

dr 


dh^ 

dr 


COS( 


dh s u 

dr 


sin (j) 


dh 


■ sm ( 


li 


dr 


COS( 


sm a. 


sin a 


( 10 ) 


The normalization is chosen so that the initial magnetic energy is Va 2 pVt^r^, where p^r { 
is unity in our scaled units. The initial field has a nonzero azimuthal component unless 


sin a = 0. Figure 1(a) shows the meridional distributions of initial angular velocity (black 


lines) and poloidal field (red lines) and figure 1(b) shows the three components of initial 
field in the equatorial plane for a = 45°. 



(a) 


(b) 


Fig. 1.— Initial conditions with magnetic obliquity a = 45°. (a) Contours of initial axisym- 
metric angular velocity (black lines) and poloidal field (red lines) in the meridional plane, 
(b) Contours of the components of initial field in the equatorial plane: B r , Bg and B $ (left 
to right). In both (a) and (b), solid lines denote positive and dashed negative. 


With the boundary conditions, the total angular momentum should be conserved, and 
with our initial conditions, it is purely axial: 

c+i 


L z = r sin 9 u^dV = (r a - r*) V] 4io / r2T k 
Jv 6 J -1 


[x)dx. 


( 11 ) 
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The radial integral can be found analytically. In our numerical calculations, e c kUj is monitored 
at each time step, and the total angular momentum is found to be conserved to one part 
in 10 8 . For the initial rotation profile (JTJ), L, = 0.9574 in our units, and the initial kinetic 
energy E 0 = 0.3191. Because of viscosity, uniform (solid-body) rotation must eventually 
be established. The final angular velocity and kinetic energy are Qqo = L z /I = 0.5714 
and Eoo = E\I2I = 0.2736, where / ~ 1.676 is the moment of inertia of the fluid shell. 
Therefore, the excess kinetic energy available for dissipation is E$ — E^ & 0.0455 ~ Eq/7. 
In §4, we compare the times required for the magnetic and excess kinetic energies to decay 
to 10% of their initial values, as functions of Rm = Re and Viz. 


3. Kinematic problem 


Before addressing the fully nonlinear problem, we study a simplified kinematic one in 
which the flow is fixed in its initial form (JTJ), and only the magnetic induction equation (J2D 
is solved. This is a linear problem and the spherical harmonics decouple. Therefore, the 
axisymmetric and non-axisymmetric fields evolve independently. 


Ra dle r (119861 )’s analysis applies here, and we summarize it. Differential rotation 


shears the axisymmetric meridional field into axisymmetric azimuthal field—the so-called 
<jj effect—with the same dependence on (r, 6) and therefore a constant characteristic decay 
time td = rl/r]. However, if diffusion is neglected, the non-axisymmetric field reverses on 
progressively finer lengthscales as time goes on. Radler demonstrates this with a cartoon 
of a field line winding up in a plane, but the point is important enough to us that we give 
a more careful argument. In ideal MHD, advection by the velocity field il(r, 9)r sin# 
preserves the meridional components B p = ( B r ,Bg ) along the flow, as can be seen by 
considering that any closed fluid contour drawn on a sphere r =constant remains on the 
sphere and encloses constant area and constant flux as it is advected; and similarly for any 
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contour on a cone 6 ^constant. Therefore, if Bp\r,9,(j>) is the meridional field at time t, 
then the evolution of this field and its derivatives is 


^SW(r, 9, 0) = Wr Bf\r, 9,0 - !J(r, 0)f] + t ^B<°> [r, 9, 0 - Si(r, 9)(], 

f) an f) 

QgBf(r, 9,0) = [r, 9,0 - n(r, 9)i] + t 9,0 - Sl(r, 9)i], 

from which it can be seen that the meridional derivatives d r B p and dgB p increase linearly 
with time unless dBp^ / dcf) = 0 (axisymmetry) or VO = 0 (solid-body rotation). Yet 
the first line above implies that J Bp ' 1 2 dV is constant. Therefore, the non-axisymmetric 
part of Bp must reverse direction on progressively smaller scales, as was to be shown. 
Finally, since the meridional field is the source of the growing azimuthal field via the term 
oc B p • VO in the induction equation, the non-axisymmetric part of must also develop 
such reversals. 


Given a small diffusivity r), it follows that the resistive timescale of the m th azimuthal 
harmonic decreases as ~ r]~ l m~ 2 | VO|~ 2 1~ 2 at late times (m 7 ^ 0) where the length 
scale is taken to be (m|VO|£) , i.e. the length scale of shear (iVOIt )" 1 divided by number 

of polarity reversals m. Once td,m t, ah components of the m th harmonic quickly decay. 
Thus the toroidal field and magnetic energy will peak at a time 

|VO |~ 2 ) 1//3 = m _ 2 // 3 (i?m / ) 1 / 3 (AO ) _1 m/0, 

( 12 ) 

Rrri (AO ) -1 m = 0 , 

where AO ~ r 0 |VO| rms is a measure of the differential rotation, and Rm! = r^Afl/r]. The 
peak magnetic energy density therefore scales as 


v mag, peak 


| (Rrri) 2 / 3 E mag:init m 4 / 3 m/0, 
I (Rm!) 2 E mag 

, init 


r\j 


rri — 0 , 


(13) 
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unless the m — 0 component of the field is isorotational. 

Figure [2] shows the time evolution of total, axisymmetric and non-axisymmetric 
magnetic energy at different Rm with Va = 0.1 and a = 45°. Initially, the axisymmetric 
azimuthal magnetic field grows linearly with time and its energy quadratically, but 
eventually all components of the held decay resistively. The non-axisymmetric energy 
grows faster and decays much faster than the axisymmetric energy. The peak value of 
non-axisymmetric magnetic energy and the time at which it is achieved are given in Table [Q 
These results can be fit by t peak oc f?m°’ 323±0 ' 004 and E magtPeak oc Am 0 753±0 ' 035 , in rough 
agreement with equations (fT2j) and (1121) . The exponent of the latter fit is influenced by 
behavior at small Rm\ between the two highest-i?m points, the logarithmic slope is 0.696. 
A toy model that we will not go into here suggests that both exponents can be expected to 
differ from their asymptotic values as Rm —>■ oo by corrections of relative size ~ Rm~ 1 ^. 
This explains the deviations of -E magiPeak from its expected asymptotic scaling rather well, 
leaving as a mystery why the t peak scaling does not deviate more than observed. 


4. Selfconsistent MHD flow 

We now study the fully nonlinear MHD flow by numerically solving both (JT|) and (j2l) . 
The importance of the magnetic forces can be estimated by comparing the peak magnetic 
energy (1T21) predicted kinematically to the energy available in the shear flow: this ratio is 


Rm 

500 

1000 

2000 

5000 

10000 

20000 

^peak 

6.5 

8.2 

10.3 

13.8 

17.2 

21.5 

-^mag, peak 

2.431E-2 

4.396E-2 

7.591E-2 

1.495E-1 

2.447E-1 

3.963E-1 


Table 1: Time to reach the peak of non-axisymmetric magnetic energy and the peak value 
of that energy for the kinematic calculations in Fig. [2] 
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Fig. 2.— Time evolution of total (solid), axisymmetric (dashed) and non-axisymmetric 
(dash-dot) magnetic energy in the kinematic problem at Va = 0.1 and a = 45°. Curves 
from bottom to top (colored black, red, blue, green, cyan and magenta) denote respectively 
Rm = 500, 1000, 2000, 5000, 10000 and 20000. 

~ Va 2 Rm 2 / 3 for m — 1. Thus if Va > RmV 1 ^ 3 , the non-axisymmetric held can be expected 


to act as a brake on the large-scale shear before it is annihilated by diffusion flSpruit 


19991 ). 


For Re -1 < Va< Rm -ly/3 , the held will be symmetrized but will eventually drive the 
how toward isorotation (in the sense of Ferraro’s law) before viscosity enforces solid-body 
rotation. Finally, if Va <C Re~ 1 , the kinematic approximation should describe the how 
reasonably well at all times. All this presumes that the axisymmetric and non-axisymmetric 
held strengths are initially comparable. The effect of varying the misalignment angle a is 
described at the end of this section. Unfortunately, numerical considerations dictate that we 
adopt Re = Rm rather than Re Rm , so we must use some care to distinguish magnetic 

from viscous effects. 


The influence of the non-axisymmetric held can also be described as a form of phase 
mixing, but one that involves azimuthal advection by the velocity held as well as Alfvenic 
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propagation along the lines. As in the axisymmetric case, fluid elements on the same 
field line can exchange angular momentum through magnetic tension, and differing rates 
of propagation on neighboring lines leads to phase mixing and damping of the associated 
Alfven waves. But since the field is non-axisymmetric, the projections of magnetic field 
lines onto the meridional may intersect. Insofar as the velocity field remains approximately 
axisymmetric and predominantly rotational, the differential rotation may bring points 
on two such lines arbitrarily close together, even though these points are initally widely 
separated in azimuth. The combination of these Alfvenic and advective processes may drive 
the entire volume to a uniform angular velocity, provided that the non-axisymmetric field 
persists. 


First we study the magnetic back reaction by varying the initial field strength. 
Figure [3] compares the time evolution of magnetic energy in the kinematic problem with 
that of the selfconsistent MHD flow at different Va. Evidently, the latter approaches 
the kinematic behavior at low Va. In the first panel, the magnetic field is so strong 
(Va — 0.1 > RmT 1 / 3 ~ 0.06) that it modifies the rotation profile significantly, as can be 


seen in Figure 4(a) The peak magnetic energy of the selfconsistent solutions is limited by 
the excess kinetic energy available in the initial differential rotation (AEq & 0.0455, SJ2J). 

Fig. H] shows the flow and field for the same simulation as in the first panel of Fig. [3] 
at the time when the non-axisymmetric magnetic energy peaks, t = 7.48. The first panel 
shows a large change in the pattern of differential rotation compared to the initial state 


(Fig. 1(a)). Meridional circulation is induced, as shown in the right panel. The winding up 


of the field is evident in the plan views shown in the last three panels. 


Next we discuss the applicability of the isorotation law B p ■ = 0. Non-axisymmetric 

field alters the conditions for isorotation. Suppose that u = RVLe^ and B = BrCr + B z e z in 
cylindrical coordinates (R, 0, z), and that Rm —> oo so that magnetic diffusion is ineffective. 
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Va=0.1 Va=0.05 




Fig. 3.— Time evolution of total (solid), axisymmetric (dashed) and non-axisymmetric 
(dash-dot) magnetic energy in both kinematic problem (black) and selfconsistent flow (red). 
Re = Rm = 5000 and a = 45°. 

In axisymmetry, the induction term V x (u x B) = R(B P ■ Vf2) e^. If this vanishes, 
then dB/dt = 0, so that isorotation can be maintained indefinitely, with different angular 
velocities on different magnetic surfaces (surfaces parallel to the held lines). When the held 
depends nontrivially on the azimuthal coordinate (0), however, the induction term becomes 

( or? or? \ 

e R + -gj-eA . 

Even if B p • VfHs initially zero, the second term on the right side will cause B p to evolve, 
so that additional constraints must be satished by the initial state to maintain isorotation. 
It is likely that there are very few if any nonaxisymmetric isorotational states other than 
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(a) 





(b) 


Fig. 4.— Sequel to Fig. Q] at the time when the non-axisymmetric magnetic energy peaks 
(t = 7.48); Re = Rm = 5000, Va — 0.1 and a = 45°. (a) Contours of axisymmetric angular 
velocity (left panel) and meridional circulation (right panel) in the meridional plane, (b) 
Contours of B r , Bg and B$ in the equatorial plane. In both (a) and (b), solid lines denote 
positive and dashed negative. 


solid-body rotation flMestel fe Weiss 


19871 ). Since the non-axisymmetric field decays much 


faster than the axisymmetric held, however, it may be that isorotation can be approached 
on intermediate timescales when the held is substantially axisymmetric but viscosity has 
not yet eliminated all shear in the velocities. Such intermediate states are possible in 
the regime Re^ 1 C VoC RmT 1 ^ 3 . Figure [5] shows the contours of axisymmetric angular 
velocity and poloidal held lines at Va = 0.1 and 0.01: the former is outside the regime 
in question, while the latter is within it. In each subhgure the four panels are taken at 
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intermediate times after the peak of the non-axisymmetric magnetic energy (see upper-left 


and bottom-right panels in figure E]) but much earlier than 12 1 Re. Figure 5(a) suggests that 


the isorotational state is destroyed by the stronger field, while Fig. 5(b) suggests that it is 
nearly achieved for the weaker held in the deeper interior. In the latter case, the angular 
velocity is nearly a function of cylindrical radius alone, as required for steady unstratified 
how when magnetic and other non-potential forces are unimportant. (Some dependence on 
z as well as R might have resulted if we had allowed for stratihcation.) The magnetic held 
lines are approximately parallel to the rotation axis at depth but not near the surface. They 
resemble the shape of the slowest-decaying resistive magnetic eigenfunction in a uniformly 
rotating and uniformly conducting shell, especially near the equator (Fig |l(a)j ). 

Thirdly we co mpare t h e rela tive strength of differential rotation and misaligned held. 


As pointed out by 


Mesteli (11999 


, §9.3), the degree of misalignment is crucial to the 


dynamics. To quantify the relative strength we compare two characteristic times. The hrst 
is the time when the excess kinetic energy drops to 10% of its initial value. As noted in 
§2, the kinetic energy decays from 0.3191 initially to 0.2736 in a state of uniform rotation 
at the same total angular momentum; thus the reduction of the difference of these by 90% 
corresponds to kinetic energy = 0.2781. The second characteristic time is that when the 
non-axisymmetric magnetic energy drops to 10% of its peak (not initial) value. Table [2] 
lists these two times for several Re = Rm and Va. On each column higher Va corresponds 
to faster decay of kinetic energy and slower decay of non-axisymmetric magnetic energy. 

A stronger held tends to convert more kinetic energy to magnetic energy, as already 
seen in Fig. 0 and the surrounding discussion above. More importantly, in the regime 
Va < 0.03 the time for the 90% drop of kinetic energy is later than the time for the 90% 
drop of non-axisymmetric magnetic energy, whereas in the regime Va > 0.05 the situation 
is reversed. Therefore, there exists a boundary in the ( Rm , Va) plane across which the 
dynamics is qualitatively changed. Below this boundary the magnetic held is symmetrized 
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(b) 


Fig. 5.— Contours of axisymmetric angular velocity (black lines) and magnetic field (red 
lines) in the meridional plane. Re = Rm = 5000 and a = 45°. (a) Va = 0.1, (b) Va = 0.01. 
In both (a) and (b) the four panels from left to right are at time = 30, 40, 50 and 100. 


before differential rotation completely decays; above it, differential rotation is suppressed 


before the magnetic held is symmet rized (if it is 


is consistent with the arguments of 


Mestell ( 19991) and 


symm etrized at a l 


Spruit 


(119991) . 


before it decays). This 


To test these arguments more quantitatively, we have interpolated in Va along each 
column of Table [2] to find the value Va CT i t at which the two times in question are equal. 
The bottom row of Table [2] indicates that Va cr it is not a monotonic function of Rm. A 
power-law fit to the final two colums (the highest Rm) yields Va CI ; t oc Am -0 ' 20 , a somewhat 
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weaker dependence than the scaling Rm predicted by the quasi-kinematic reasoning at 
the beginning of this section. We return to this point in §3 below. 

To end this section, we briefly discuss the effect of varying the angle a between the 
rotational and magnetic axes. Since, as we have seen, there are strong differences between 
the axisymmetric cases (where a = 0) and those for which a = 45°, it is reasonable to 
expect that the strength of the interaction between the flow and the field should increase 
continuously with this angle up to a = 90°. This expectation is tested in Figure EJ which 
shows the evolution of the magnetic energy and Ohmic dissipation for several values of 
a. A larger angle leads to higher energy and dissipation, as expected. In fact, for the 
cases shown in Fig. [6] and tabulated in Table |3j the peak of the non-axisymmetric energy 
increases linearly with sin 2 a (correlation coefficient 0.9999). This would not be surprising 
in the kinematic problem, since the initial amplitude of the non-axisymmetric component is 
proportional to sin a, but the calculations shown in Fig. [6] include magnetic backreaction. 
Of course the energy E mSLg!peak (and somewhat more, because dissipation has already begun 
to act at the time of the peak) must come at the expense of the differential rotation. 


5. Discussion 

In this work we have explored numerically the interaction of differential rotation and a 
misaligned magnetic field. In the kinematic limit we verify the Rm 1 ^ 3 law for the time and 
amplitude at which the non-axisymmetric azimuthal field reaches its peak and then decays, 
whereas the corresponding scaling for the axisymmetric field is O(Rm). In selfconsistent 
calculations where the flow reacts to magnetic forces, a sufficiently weak magnetic field 
behaves approximately as in the kinematic problem, becoming axisymmetric before the 
rotation becomes uniform, and a nontrivial state of isorotation—one in which different 
magnetic surfaces have different angular velocities—is approximately achieved near the 
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500 

1000 

2000 

5000 

10000 

0.01 

20.6/17.9 

39.5/20.2 

75.0/23.9 

159.4/31.1 

258.8/38.2 

0.02 

19.8/18.1 

36.5/20.5 

64.1/24.2 

115.6/31.6 

185.9/38.9 

0.03 

18.4/18.4 

32.0/21.0 

50.7/25.1 

72.5/32.3 

60.4/38.7 

0.04 

16.5/18.9 

26.3/21.7 

35.3/26.0 

34.2/32.7 

30.2/45.0 

0.05 

14.4/19.5 

20.0/22.9 

21.8/27.0 

23.9/33.4 

25.4/47.0 

0.1 

8.5/21.3 

9.8/26.2 

11.8/31.0 

19.3/37.0 

— 


0.0300 

0.0461 

0.0464 

0.0414 

0.0359 


Table 2: The time when the kinetic energy drops 90% of its initial value (on the left of slash) 
and the time when the non-axisymmetric magnetic energy drops 90% of its peak value (on 
the right of slash). The top row indicates the values of Re = Rm, the left column Va and 
the bottom row Va cr i t . a = 45°. The calculation at Re = Rm = 10000 and Va = 0.1 was 
incompletely resolved. 


a 

0° 

30° 

45° 

60° 

90° 

-^mag, peak 

0 

7.019E-4 

1.387E-3 

2.057E-3 

2.711E-3 


Table 3: Peak non-axisymmetric magnetic energy versus magnetic obliquity a, following 
Fig ED Re = Rm = 5000 and Va = 0.01. 
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Fig. 6.— Time evolution of magnetic energy (left panel) and Ohmic dissipation (right panel) 
at different magnetic obliquities. Re = Rm = 5000 and Va = 0.01. Black, red, blue, green 
and cyan lines denote respectively a = 0°, 30°, 45°, 60° and 90°. Solid, dashed and dash-dot 
lines denote respectively the total, axisymmetric and non-axisymmetric energy or dissipation 
(a = 0° and a = 90° are entirely axisymmetric and non-axisymmetric, respectively). 

rotation axis. For stronger fields, there exists a boundary in the plane of dimensionless 
initial magnetic field strength Va versus Rm above which differential rotation is suppressed 
before the nonaxisymmetric magnetic field decays. This boundary is better described by 
Va cr it oc Rm~ 0 ' 20 than by the expected scaling Rm _1//3 , at least over the range explored 
by our simulations (Rm < 10 4 ). We do not understand this quantitatively for our actual 
simulations, though according to the arguments given below, one expects the dependence 
of Va CI ^ on Rm to weaken at very large Rm because of reconnection and perhaps other 
nonlinear processes that assist in the destruction of the nonaxisymmetric field. We have 
also verified that a larger misalignment angle (a) between the rotational and magnetic axes 
leads to stronger exchanges between the flow and the field, in proportion to sin 2 a (Fig. |6|). 

In our calculations the magnetic Prandtl number Pm = Rm/Re = v/r\ has been fixed 
at unity for numerical convenience. With the code and computational resources available 
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to us, fully resolved simulations are practical only if the larger of Re and Rm is < 10 4 


Therefore, we approach ideal MHD best by ta 


of ma gneti c dyn amos (e.g. 


dng Re = Rm. Direct numerical simulations 


Schekochihin et al. 


(iFromang et al. 


2007}) and magnetorotational turbulence 


20071. e.g.) often find that Pm influences large-scale properties of the 


flow, at least within the computationally accessible range of Re and Rm. Whether Pm 
remains macroscopically important at much larger Rm and Re, as in real stars and in the 


i nterst ellar medium, c an o n 


(IKulsrud & Anderson 


1992 


y be addressed at present via idealized moc 


Malvshkin & Boldyrev 


2009 


Schober et al. 


els of turbulence 


2012 ). It is worth 


noting, however, that in comparison to Re and Rm themselves, Pm is not so far from unity 
in the applications of interest to us. At r = O.5R 0 in the Sun (where T « 4 x 10 6 K and 
p ~ 1.3gcm -3 ), one estimates from the classical Spitzer formulae that Pm pz 10~ 2 . For 
white-dwarf interiors at T ~ 10 6 -10'K and p ~ 10 6 gcm~ 3 , we estimate! Pm ~ 10 2 -10 3 . 

Another simplification we have made is to neglect stratification, an important effect 
in the outer parts of white dwarfs, where there are significant entropy and composition 
gradients, and in the radiative zone of the Sun. As noted in §1, stratification would 
allow a larger range of isorotational final states in hydrostatic equilibrium. The ratio of 
Brunt-Vasala frequency to rotation frequency (N/Q) is important for the growth rates of 


magnetic instabili 


itself (Spruit 


ies that may gradually modify the differential rotation as well as the field 


1999!) . However, we do not expect that stratification has a major effect on the 


interaction between a non-axisymmetric field and differential rotation (e.g., the dependence 
of Va CT it on Rm), at least on the dynamical timescales simulated here. Nor have we allowed 


using the results of Yakovlev & Urpin fll980h as encoded in the FORTRAN programs 


at www.ioffe.ru/astro/conduct, which return electrical and thermal conductivities. We 
roughly translate the latter to a kinematic viscosity by multiplying by rrii/k B p, rrii being the 
ion mass, presuming that the material has not crystallized. 
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for compressibility or density gradients. All of these could be the subject of future work, 
perhaps based on realistic models of representative stellar models. 

Notwithstanding these caveats, we tentatively suggest some implications for real 
stars. As discussed in §1, this work has been motivated in part by tidal interactions 
between white dwarfs in binary orbits decaying under the influence of gravitational 
radiation. The calculations we have performed, however, are not directly applicable to that 
astrophysical problem. For one thing, the magnetic R eynolds n um ber o f a white dwarf 


198.0), radius 


is enormous: Rm ~ 10 1 ' for diffusivity r] cs 1 cm 2 s _1 (lYakovlev fc IJrpin 
R ~ 0.01 Rq « 7 x 10 8 cm, and angular velocity = 27r/lmin. Solar flares and other 
more theoretical evidence suggest that under such nearly ideal conditions, the destruction 
of magnetic energy proceeds more rapidly than by linear resistive diffusion. It is plausible 
that nonlinear dissipative processes are already important at the modest values of Re and 
Rm ~ 10 3 in our simulations. This may account for the difference between the predicted 
and observed scalings of Va CTlt with Rm. While a variety of collisionless-plasma effects have 
been invoked to explain the observed rate magnetic energy release in solar flares, recent 
theoretical work indicates that even in classical resistive MHD, oppositely-directed field 
lines may approach one another and annihilate at a speed V Tec ~ /Va that is independent of 
the true microscopic diffusivity ( r /), provided that the Lundquist number S = LVa/t) > 10 4 , 
where Va is the physical Alfven speed and L is an appropriate macroscopic lengthscale, such 
as the length of the recon necting current sheet; the present estimate for the dimensionless 


factor is / ~ 0.02 (iLoureiro et ah 


2012 


, and references therein). Mechanisms by which 


Vrec/Vx becomes independent of r] (or perhaps depends upon it only very weakly) are said 
to provide “fast reconnection”. 


In our problem, Va should be based on the tightly wound azimuthal field, which scales 
with time t as Va ~ VaRQAQt as explained in §|3l while the lengthscale along field lines is 
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L ~ ro- Thus S ~ VaRm!. The local timescale on which the field is destroyed now becomes 
l/V rec rather than Z 2 / 77 , with l ~ r 0 /(tAf2) being the cross-field lengthscale on which the 
field reverses. Equating this local timescale to the winding time t leads to a revised estimate 
of the time at which the magnetic energy should reach its peak: Af2t p k ~ (/ Vafi/Afi ) -1 / 3 
rather than RrrT 1//3 as before. Defining Va cr i t so that the magnetic energy at the peak 
is equal to the excess kinetic energy A E in the initial differential rotation leads to 
Ta cr it (D//AD) 1 / 2 (AE/5-E' 00 ) 2 , where E 00 = L 2 /2J is the final rotational energy, and the 
factor 1/5 arises from I/Mr\ cs 2/5 as for a full sphere of mass M and constant density. 
Finally, setting / ~ 0.02, AD/f2 ~ 1, and A E/E^ ~ 1/6 yields Va CTit ~ 0.011. In short, if 
we were able to extend our calculations to Re = Rm —> 00 with the same initial rotation 
profile and magnetic geometry, then we would expect Va cr i t to asymptote to ~ 10 -2 because 
of fast reconnection. It is unlikely that our simulations achieve fast reconnection, however, 
because the Lundquist number at the peak is 


5 


pk 


Rm 


—sr) 


where we have put Va = 0.04 and Af2/f2 = 1 in the final estimate. Thus our largest-ZZm 
simulations probably imdergo stable Sweet-Parker reconnection, for which V rec oc V A / and 
an argument along the lines above predicts Va cr ; t oc RrrT 1 ^^. I 11 astrophysical applications 
where S 10 4 , it is likely that fast reconnection dominates so that Va cr it becomes 
independent of Rm. 


In an effort to cast dou bt on conte mpor ary h elioseismological evidence for differential 


rotation in the solar core, 


Meste l fe Weiss ( 19871 ) estimated that a meridional field 


B P > 0.03G would be sufficient to establish and maintain isorotation in the radiative zone. 
They admitted that isorotation does not require uniform rotation but went on to speculate 
that the latter would result if the magnetic field were “even slightly non-axisymmetric.” 
We are perhaps in a position to quantify this. From the results of Fig. [ 6 ] and surrounding 
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discussion, Mestel and Weiss’s estimate (of conditions sufficient to establish uniform 
rotation) can be sharpened to B P sin a > 0.03G, where a is the magnetic obliquity. 

But this presumes that the differential rotation does not symmetrize the field (reduce a 
to zero) before solid-body rotation is established. If the field in the core is frequently 
regenerated by a dynamo process, we know of no way to constrain its obliquity. But if 
it is a fossil field established in the early history of the Sun, then the free-decay problem 
we have studied may apply. The magnetic Reynolds number of the radiative zone is 
Rm & = Rlfl & /r] ~ 10 14 , where we take R c ~ 0.7R© for the outer radius of this zone 
and evaluate 77 for T = 4 x 10 6 K and p = 1.3gcm -3 (i.e., conditions at O.5R 0 ). If 
kflcrit = RmR 1 ^ 3 , then to avoid symmetrization, the field and obliquity would have to satisfy 
B P > 10G, considerably stronger than the estimate above but still quite modest. However, 


the Lundquist number R c fl Q /r] of the radiative zone is on the order of 10 8 (R/1G), so 
according to the discussion above, we expect to be in the regime of fast reconnection where 
Vocrft becomes independent of Rm. Adjusting mom ent o f inertia in the argument above for 
the actual density profile of the solar core ( Bahcall et al. 199hl ). we estimate Va cr ; t ~ 0.004, 


which corresponds to B P sin a ~ 2kG. It is still debated whether a fossil magnetic field in 
the solar radiative core can enforce solid-body rotation, however. A crucial issue is whether 
such a field can remain closed within the core, because if it were to conn ect to t he convectio n 


zone, then by Ferraro’s L aw t he c ore a. 


1998 


MacGregor fe Charbouuean 


Acevedo-Arreguin et al 


1999 : 


so shoul d rotat e differentially (Gough & 


Brun & Zalin 


2006 


Garaud fe Gar and 


Mclntvre 


20081 : 


2013). 


A second difference between our simulations and the binary-white-dwarf problem is 
that we have studied the free decay of a pre-existing profile of differential rotation rather 
than gradual acceleration by a tidal torque. (This allowed the outcome to be discerned with 
less CPU time.) In the tidal-binary scenario, the system starts with an orbital period of a 
few hours, which is short enough to lead to merging of a pair of 0.7 M© white dwarfs within 
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10 9 yr (31]). Assuming that the tidal torque is absorbed in the outer layers of the star, the 
magnetic stress necessary to maintain synchronous rotation of the interior is 


Br-Brl 


i dn 
B?~dt 


P 


- 11/3 


MG 2 , 


1 miry 

in which I fa 0.2 MR 2 is the moment of inertia of a cool white-dwarf model at this mass. 
Thus for example, if the meridional field is ~ 10 kG, a plausible upper limit for most 
white dwarfs, then only a slight bending of the lines is needed to maintain synchronism 
when the period is an hour or more. When the period is about a minute, however, 

B$ r^i 10 4 H r 100 MG would be required since differential rotation alone will not increase 
B r . While this is not outside the range of surface fields observed in some white dwarfs, a 
non-axisymmetric field so tightly wound would have a resistive time ~ 10~ 8 times smaller 
than that of the most slowly decaying magnetic eigenmode (~ 3 x 10 9 yr); and magnetic 
reconnection would probably act on even shorter timescales, as discussed above. Hence 
an axisymmetric field might be expected, unless perhaps dissipation in the interior leads 


to nonaxisymmetric turbulence that produces a magnetic dynamo. In fact, as 


Spruit 


(120021 ) has emphasized, a predominantly azimuthal axisymmetric field should be subject 


to nonaxisymmetric Tayler instabilities, though in strongly stably-stratified regions these 
instabilities will be concentrated toward the poles, whereas it is near the equator that 
meridional field is most needed to transmit tidal torq ues. Al so, while some numerical 


evidence for such a dynamo has been claimed by 


Brai thwaite ( 2006. 1. it is stil 


whether there exists dynamo action in stably stratified zones (e.g., 


Zahn et al. 


debated 


2003). If 


it exists, such a dynamo might maintain B r in a constant ratio to B allowing a lower 
overall stress for the same torque. By our estimates, a radial field of order 1 MG would also 
begin to have a significant effect on the dispersion relation of g-modes at t he compositional 


interface where they are tidally excited in the models of 


Ful le r fe Lail (2012). 


In summary, while there is little doubt that the interactions among tidal torques, 
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differential rotation, and magnetic fields are important for coalescing white-dwarf binaries, 
we are far from being able to predict the tidal heating and luminosity of such systems in 
the last stages of inspiral. 

We thank Prof. Rainer Hollerbach for help with the calculation of energy and 
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materially improved the presentation. JG thanks the Australian National University and 
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